R/Function 9 oceanVulnerability v.0.04.r

#' Function 9 oceanVulnerab: link hotspots of biodiversity to projected climatic impacts
#'
#' @description access and links climatic projections to hotspots of biodiversity to unravel potential impacts of climate change 
#' @param biodiversity_grid raster layer with any biodiversity metric previously computed with oceanDiversity
#' @param biodiversity_metric the name of the biodiversity metric
#' @param reshape_climatic_layer whether climatic layer resolution must be changed to match biodiversity grid resolution
#' @param new_cell_size define desired cell size of the climatic layer for the analysis
#' @param map_climtatic_impacts whether we want the resulting spatial grid of future trends to be mapped
#' @param plot_histograms whether histograms of expectec temperature increases at hotspots are to be mapped
#' @return Object of hotspots coordinates and climatic projections according to IPCC scenarios RCP 2.5, RCP 4.5, and RCP 8.5
#' @details This function allows accessing and exploring expected changes environmental parameters 
#' @export


oceanVulnerability = function (biodiversity_grid,
                          climatic_layer =NULL,
                          get_climatic_layers = T,
                          IPCC_scenario = "RCP85",
                          reshape_climatic_layer = F,
                          new_cell_size = 1,
                          plot_histograms = T,
                          map_climtatic_impacts = T,
                          low_color="steelblue", 
                          mid_color="gold",
                          high_color= "firebrick", 
                          col_steps=20) {
  

if(reshape_climatic_layer){
    
    factor_conver = new_cell_size / res(climatic_layer)[1]
    climatic_layer <- aggregate(climatic_layer, fact=factor_conver)
    
}
  
  hotspots = oceanHotspots(biodiversity_grid, only_hotspots = T, map_hotspots = F)
  
  hotspots = as(as(hotspots, "SpatialPointsDataFrame"), "data.frame") 
  
  colnames(hotspots) = c("hotspot", "x", "y")
  
  xy = SpatialPoints(hotspots[,c("x","y")])
  projection(xy) <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
  
if(get_climatic_layers){
    
    climatic_layerRCP = oceanFuture(IPCC_scenario = IPCC_scenario)
    
    hotspots$future_SST = extract(climatic_layerRCP, xy)
}
  
color.gradient <- function(x, colors=c(low_color,mid_color,high_color), colsteps=col_steps) {
return(colorRampPalette(colors) (colsteps) [ findInterval(x, seq(min(x),max(x), length.out=colsteps)) ] )}
  
  
if(plot_histograms){
      
      par(mfrow=c(1,2)) 
      
      hist(hotspots$future_SST, breaks=15, main=paste(IPCC_scenario,"IPCC"),
                 col="red",xlab="Expected Increases in Tª(ºC)", 
                 ylab="Hotspots of Biodiversity")
      box()
      
    }
    
  
  
if(map_climtatic_impacts){
    

palette(color.gradient(1:20))
    
plot(biodiversity_grid, col="white", legend=F, main=paste(IPCC_scenario,"IPCC"))

maps::map("world",add=T, fill = T,bg="grey20",col="grey30")
points(hotspots[,c("x","y")], pch=21,cex=1.5, 
           bg=factor(hotspots$future_SST))
    
  }

par(mfrow=c(1,1)) 

  return(hotspots)
  
}
monteroserra/rOceans documentation built on Dec. 23, 2024, 2:47 p.m.